Convergence of density-matrix expansions for nuclear interactions 
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We extend density-matrix expansions in nuclei to higher orders in derivatives of densities and test 
their convergence properties. The expansions allow for converting the interaction energies charac- 
teristic to finite- and short-range nuclear effective forces into quasi-local density functionals. We 
also propose a new type of expansion that has excellent convergence properties when benchmarked 
against the binding energies obtained for the Gogny interaction. 
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The description of finite many-body systems is a con- 
tinuous exciting challenge in physics. In view of the ex- 
ploding dimensionality of many-particle Hilbert spaces, 
direct approaches are however often untractable. Instead, 
density functional theory (DFT) provides an alterna- 
tive by showing that a much simpler description in terms 
of the density alone is in principle possible. The main 
challenge is finding or deriving the best density function- 
als that are able to consistently describe large classes 
of physical systems. In this work we address a related 
question, namely to what extent the physics of the exact 
nonlocal exchange potential can be captured by an ap- 
proximate density functional? This question is not only 
of fundamental importance but also of great practical in- 
terest, and attempts to solve it date back to the early 
1950s, involving the Slater as well as the optimized effec- 
tive potential approximations [H. 

For nuclear-physics applications, Negele and Vautherin 
developed an alternative method known as the density- 
matrix expansion (DME) [5], which takes advantage of 
the short-range of the nuclear two- and three-body forces. 
It is particularly convenient as it was originally con- 
structed to generate quasi-local density functionals of 
a form similar to the successful phenomenological Skyrme 
functionals used extensively today Ijl . The method have 
recently gone through a revival [5"9] and appears to be 
a promising starting point for the construction of ah ini- 
tio nuclear DFT models 0, S] . The idea is to start from 
effective forces, derived from first principles, and then 
transform the expression for the interaction energies into 
quasi-local density functionals using the DME. Because 
the starting point is an effective force, this method is able 
to take correlations into account in a way that goes much 
beyond the Hartree-Fock treatment of the bare interac- 
tion. 

In conjunction with these ab initio methods, there are 
empirical efforts to improve nuclear density functionals, 
such as the systematic method we recently proposed to 
construct nuclear densit y fu nctionals using higher-order 
derivatives of densities jlO| . In such an approach, the 
DME is helpful as it can serve as a method of classifying 
different terms according to their importance, similarly 



to what is available in the effective-field-theory power- 
counting scheme [ll|. 

In this Letter, we explore the ideas of the DME to 
transform the direct and exchange energies resulting from 
the finite-range Gogny interaction into quasi-local den- 
sity functionals. We then extend these methods beyond 
the second order expansions used so far to make the first 
tests of convergence and accuracy in function of the ex- 
pansion order. 

We begin our analysis by recalling that the nuclear 
one-body density matrix p (tictiTi, r2cr2T2), which de- 
pends on the space-spin-isospin nucleon coordinates rur, 
can be separated into four nonlocal spin-isospin densities 
/9*(ri,r2) as 



p(riCTlTl,r2Cr2T2) = 

V + t r 
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[<,^,p*(ri,r2)]°l , (1) 



where we sum over the spin-rank v = Q and 1 (scalar and 
vector) and isospin-rank t = Q and 1 (isoscalar and isovec- 
tor) spherical tensors of the Pauli matrices, Ov and 
r*, respectively. We denote vector (isovector) coupling 
by the standard square brackets with subscripts (super- 
scripts) giving the values of the total spin (isospin) . 

For a local interaction, the energy corresponding to 
the direct (Hartree) term depends on the local densities 
p* (r) = p* (r, r) in Eq. ^ as 
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[p*(ri),p*(r2)]" ,(2) 



where V^" (»"i,»"2) denote the local direct two-body po- 
tentials. As discussed in Ref. Q, for short-range inter- 
actions and sufficiently smooth local densities one can 
employ simple Taylor expansions of densities to obtain 
the energy density in the form of a gradient expansion. 
As in Ref. [13], we use the spherical tensor representa- 
tion of derivatives, which allows for the most economical 
classification of terms in the expansions. 

Introducing the standard average R— ^ (ri + r2) and 
relative r = (ri — r2) coordinates, we have the following 
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TABLE I: The direct end exchange energies in ^"^Pb calcu- 
lated with the finite-range Gogny DIS interaction 115] in the 
scalar {v — 0), vector {v — 1), isoscalar (t = 0), and isovector 
{t = 1) channels. All energies are in MeV. 
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expansion of the local densities, 

71 

nL 

where, for clarity, the spin and isospin projection quan- 
tum numbers, ruv and mt, are not shown, V denotes 
the gradient operator with respect to variable R, DnL 
denote the spherical-tensor derivatives of order n and 
rank L Yl if) are the standard spherical harmon- 
ics, and a„L are simple numerical constants resulting 
from the recoupling of the Cartesian nth order deriva- 
tives (^^J" • into the spherical-tensor forms [isj . 

Expanded in this way, and after a partial integration 
and recoupling needed to make all derivative operators 
act on the same density, the direct interaction energy of 
Eq. ([2]) can be written as a gradient expansion. 
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In this final form, the coupling constants C^j^ are num- 
bers obtained as moments of the direct potentials up to 
the order of the expansion N. This order is defined by 
keeping all terms in the product that fulfill the condition 
n + n' < N , where n and n' denote the orders in the 
expansions of the two densities [compare with Eq. ([3])]. 

To test the quality and convergence properties of this 
expansion, we generated a set of realistic densities for 
doubly magic nuclei by solving the self-consistent equa- 



tions for the nuclear SLy4 functional For these 

densities, the sums in Eq. ^ were evaluated to differ- 
ent orders by calculating the higher-order derivatives of 
densities needed using the program hosphe (vl.OO) [l^ . 
The coupling constants in Eq. ^ were calculated from 
moments of the finite-range part of the Gogny DIS in- 
teraction Then for the same set of densities, exact 
energies were determined using the capability of treating 
Gaussian interactions with the code hfodd (v2.45h) [lij . 



As an example, in Table |T] we show the exact energies 

calculated in ^"^Pb for the four spin-isospin channels. In 

the ground states of even-even nuclei, the vector terms 
IOFt; 
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FIG. 1: (Color online) Precision of the gradient expansions 
(Q of the direct interaction energies calculated for the finite- 
range Gogny DIS interaction [1^. The nine nuclei used for 
the test are *He, i^^O, «'*«Ca, ^^.^SNi, ^"0,132^^^ 208p^_ 
Positive (fourth order) and negative (second and sixth orders) 
errors are shown with open and closed symbols, respectively. 



of the direct energies vanish due to the conserved time- 
reversal symmetry, and thus our calculations here only 
test the scalar terms of the expansion in Eq. In 
Fig. [T] we show percentage deviations between the gradi- 
ent approximations and the exact results. It is extremely 
gratifying to see that in light nuclei the sixth-order ex- 
pansion allows for reaching the precision of 0.1-0.01% or 
0.1 MeV. Even more important, in each higher order the 
precision increases by a large factor, which is character- 
istic to a rapid power-law convergence. One also sees 
that in heavier, isospin-unsaturated nuclei the precision 
slightly deteriorates, probably owing to a slower conver- 
gence for the isovector densities, which fluctuate more 
than the isoscalar ones. In comparison, models based on 
the best tuned Gogny forces give RMS deviations of 0.798 
MeV [ftI for the description of experimental masses and 
to go much beyond this precision in the expansion is thus 
not so useful. 

One can, of course, quite easily calculate the direct 
terms exactly, but still it is interesting to see that a few 
of the lowest-order coupling constants appearing in the 
expansion are able to capture the essential physics con- 
tained in the nuclear finite-range force. For the exchange 
(Fock) terms whereto we turn now, the situation is signif- 
icantly different, because here the expansion into quasi- 
local densities dramatically increases the simplicity and 
efficiency of calculations. 

In order to deal with the exchange term, we propose a 
new approximation, which we call a damped Taylor (DT) 
expansion, and we define it in the following way: 
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a= (p*(ri,r2) -p*(ri,r2)) 
I 



= Y.<Lir)[YLif),pUAR)]o- (5) 

ri=r2 = R nLm 



The Taylor expansion of the term in the square brackets 

is obtained by using e*'' '^ = X^m 7^ lir-kj , where 

k — (Vi — V2) is the standard relative-momentum 
operator. The quasi-local densities, are defined as 



pliLviR) = KnLpl{ri,r2) 



(6) 



where the higher-order spherical-tensor derivative oper- 
ators KnL UM built from the first-order operators 
k. Because of the simple form of the expansion, explicit 
formulas for functions tt^ (r) can be derived and will be 
presented elsewhere [l^l- 
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FIG. 2: (Color online) Comparison of a Taylor expansion and 
the DT expansion (with = and a = A/kp), applied to 
the SNM nonlocal density in second, fourth, and sixth orders. 
The figure is drawn for k% = 1.35 fm^^ corresponding to 
the equihbrium density of SNM d, The shaded regions 
indicate the integration ranges needed to obtain 1, 0.1, 0.01 
and 0.001 % accuracy for the SNM exchange energy with the 
Gogny DIS interaction. 

The DT expansion has the potential to be more accu- 



rate than the original Negele-Vautherin (NV) expansion 
because it avoids the approximation of angular aver- 
aging the density matrix in the nonlocal r direction. It is 
generated around a model density p* and is also damped 
by a Gaussian factor with width a. In the limit of van- 
ishing model density and infinite damping distance, the 
DT becomes an ordinary Taylor expansion. This Taylor 
expansion applied to the uniform system of symmetric 
nuclear matter (SNM) is illustrated in Fig. [5] It works 
well in a region around r — but diverges for larger rel- 
ative distances. These bad asymptotics can however be 
cured by using a finite damping width a, as also demon- 
strated in Fig. [2] By employing the condition that it 
should work equally well in SNM for all values of Fermi 
momentum kp, we estimate this width as a = A/kp. 
Assuming that the densities of real nuclei have similar 
behavior, we retain this estimate of a in the calculations 
presented below. 

The model density can be chosen in different ways by 
the requirement that the truncated expansions must be- 
come exact in certain limits. The model used in the fol- 
lowing leads to exact values in SNM, and is obtained 
by assuming that the density in each point R has the 
SNM nonlocal dependence, p*(ri,r2) = pl{R) ^^^^jf^- 
For large relative distances, when the expansion around 
r — can no longer be trusted, the bracketed term in 
Eq. ([5]) is damped to zero and the asymptotic behav- 
ior is then given entirely by this model part. In this way 
the parameter a interpolates between placing trust in the 
model part or in the expansion part, and the optimal bal- 
ance is determined by the range in which the expansion 
converges, as estimated above. 

By inserting the DT expansion of Eq. ([S]) into the ex- 
pression for the exchange energy, we obtain the quasi- 
local expansion 
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where the quasi-local densities PnLyj{R)^ which were in- 
troduced in Ref. 10], are equal to those of Eq. ([6]) with 
the angular momenta L and v coupled to the total an- 
gular momentum J. The coupling-constants Cl^^vj 
given by the integrals of the local exchange two-body 



potentials V°^'^ (ri,r2) = V^f^ (r) and functions ttJ^ (r) 
of Eq. ([5]). They depend on kp through the model 
density and damping width, which in the standard lo- 
cal density approximation leads to a density dependence 
kp = (3pO/27r2)V3 0. 
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FIG. 3; (Color online) The RMS deviations between the ex- 
act and approximate exchange energies of Eq. ((Tj, calculated 
for the nine nuclei listed in the caption of Fig. [T] The four 
panels show results obtained in the four spin-isospin channels 
labeled by Vvt- Results obtained in the present work (DT) are 
compared with those corresponding to the NV 0] and PSA 
[3| expansions. 



In Fig. [3] we show results obtained with the DT ex- 
pansion compared to the exact values, similar to those 
listed in Table HI The comparison is made by considering 
the RMS deviations between the exact and approximate 
exchange energies, and is also shown for the NV 2] and 
phase-space averaging (PSA) Q expansions. As seen in 
this figure, the NV Bessel-function expansion used in the 
scalar channels works very well, whereas the alternative 
method used for the vector channels [see Eq. (4.23) in 
Ref. [3]] does not improve when evaluated beyond the 
second order. 

We employ the recently proposed PSA method in the 
INM-DME version and generalize it to higher orders 
by incorporating the model density in the same way as 
in Eq. ([5]). Other modifications are probably also possi- 
ble and could be studied as well. At second order, our 
prescription reproduces that proposed in Ref. with 
the PSA performed over SNM. Fig. |3] shows that this is 
the best expansion at second order, where the total RMS 
deviation for the exchange energy is 14.09 MeV. When 
used for the vector part, the PSA expansion is similar to 
the NV expansion and also has problems when evaluated 
beyond second order. 

As seen in Fig. |31 the DT expansion exhibits excel- 
lent convergence properties both for the vector and scalar 
parts, with every next order improving the results by 
large factors. Note that in the vector channels, the exact 
exchange energies of Table |T] are quite small, so at second 
order their relative description is poor. In contrast, at 
higher orders their description becomes rather accurate. 
The results obtained for the DT expansion are not so 
sensitive to the value of the damping width a, and values 
in the range of a « S/kp - 5/kp all lead to improvements 
over the NV expansion. 

The total RMS deviations obtained in the DT case for 



the exchange energies at fourth and sixth orders are 1.54 
and 0.36 MeV, respectively. Calculations using a fixed 
value of kp = 1.35 fm^^ give results that are not signifi- 
cantly different [13|. They lead to substantial simplifica- 
tions, because then the density dependencies of coupling 
constants can be neglected. 

In summary, we presented the first analysis of density- 
matrix expansions extended beyond second order, with 
derivatives of densities included up to sixth order. We 
also proposed a new type of expansion that has excel- 
lent convergence properties. The results demonstrate 
that its possible to construct expansions into quasi-local 
density functionals that converge for the description of 
low-energy nuclear observables owing to the short range 
of nuclear forces. These results also provide a baseline 
for phenomenological adjustments of quasi-local higher- 
order nuclear functionals. 
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